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Flow  separation  on  the  lifting  surfaces  of  a  vehicle  at  high 
angle  of  attack  is  always  complicated  by  a  certain  degree  of 
unsteadiness/  but/  when  the  vehicle  itself  is  undergoing  unsteady 
motion  or  deformation/  or  if  it  enters  a  different  flow  field 
rapidly/  then  the  complexity  of  the  separated  flow  is  even 
greater/  and  culminates  in  the  phenomenon  of  dynamic  stall.  If 
the  angle  of  attack  oscillates  around  the  static  stall  angle/  the 
fluid  dynamic  forces  and  moments  usually  exhibit  large  amounts  of 
hysteresis  and  a  condition  of  negative  aerodynamic  damping  often 
develops  during  part  of  the  cycle.  This  can  lead  to  the  condi¬ 
tion  of  flutter  in  a  single  degree  of  freedom  oscillating  rigid 
body  motion.  (Normally/  in  attached  flow,  flutter  only  occurs 
when  the  body  motion  includes  multiple  degrees  of  freedom;  e.g., 
combined  bending  and  torsion  of  an  aircraft  wing.)  During  a 
rapid  increase  in  angle  of  attack ,  the  static  stall  angle  can  be 
greatly  exceeded/  resulting  in  excursions  in  the  dynamic  force 
and  moment  values  that  are  far  greater  than  their  static  counter¬ 
parts.  The  consequences  of  dynamic  stall  are  far-reaching  and 
lead  to  such  problems  as  wing  drop,  yaw  (sometimes  leading  to 
spin  entry)/  wing  rocking  and  buffeting  as  well  as  stall  flutter. 

Although  a  great  deal  has  been  learned  about  dynamic  stall 
characteristics — mainly  through  experimental  observation — there 
is  not  at  this  time  a  completely  satisfactory  theoretical  method 
(1),  (2)  for  predicting  the  dynamic  stall  characteristics  for  new 
untested  shapes,  even  for  the  two-dimensional  case.  Moreover, 
quantitative  comparisons  of  experimental  test  data  on  new  geome¬ 
tries  can  be  obscured  by  the  effects  of  three-dimensional  wind- 
tunnel  interactions,  wall  interference  and  experimental  uncer¬ 
tainties  (3).  In  the  present  work  a  possible  theoretical  ap¬ 
proach  is  examined  for  predicting  dynamic  stall  characteristics. 
The  approach  combines  an  unsteady  time-stepping  method  (4)  with  a 
steady  inviscid/viscid  iterative  code  (5)  that  includes  an  exten¬ 
sive  separation  model.  The  latter  has  proven  very  successful  in 
the  steady  case.  Both  codes  are  applicable  to  general  three- 
dimensional  shapes  and  have  been  developed  from  the  same  basic 
panel  method  (6). 


Extensive  investigations  of  the  dynamic  stall  characteris¬ 
tics  of  airfoils  oscillating  in  pitch  have  been  reviewed  by 
McCroskey  (1),  (2).  In  practical  aerodynamic  environments,  the 
unsteadiness  can  be  a  combination  of  several  motions.  Unsteady 
motions  other  than  pitching  have  been  investigated  by,  among 
others,  Liiva  et  al.  (7),  Lang  (8),  and  Francis  (9),  Maresca  et 
al.  (10),  and  Saxena  et  al.  (11).  In  these  references  the  phen¬ 
omena  of  plunging,  flap,  spoiler  and  rectilinear  oscillations 
were  examined.  There  are  very  few  theoretical  approaches,  how¬ 
ever,  and  in  a  recent  review,  McCroskey  (2)  concludes  that  at  the 
present  time  the  engineer  who  needs  answers  should  turn  to  one 
of  the  empirical  correlation  techniques,  even  though  these  are 
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not  completely  satisfactory  and  only  supply  broad  details  of 
force  and  moment  characteristics.  The  method  of  Ericson  and 
Reding  is  perhaps  one  of  the  most  comprehensive  of  these  methods. 
Their  latest  paper  (12)  incorporates  a  number  of  findings  from 
the  systematic  experiments  of  Carr  et  al.  (13). 

Early  theoretical  approaches  to  dynamic  stall  addressed  the 
deep  stall  aspect  which  is  dominated  by  the  passage  of  vortices 
shed  from  the  vicinity  of  the  leading  edge.  Ham  (14)  used  just  a 
thin  plane  model  of  the  airfoil.  Later  work  by  Baudu  et  al.  (15) 
extended  the  basic  model  to  thick  sections  using  a  panel  method. 
The  main  drawback  with  the  approach  is  that  crucial  assumptions 
regarding  the  location  and  time  of  vortex  shedding  have  to  be 
made  in  order  to  perform  the  calculations.  Also,  the  results 
(from  (15))  are  sensitive  to  (a)  the  angle  at  which  the  vortex 
path  leaves  the  surface/  (b)  the  time  at  which  vortex  emissions 
terminate  so  as  to  start  the  reattachment  process  and  (c)  the 
viscous  diffusion  of  the  free  vortices. 

Calculations  of  the  characteristics  of  the  thin  boundary 
layer  in  the  attached  flow  regions  of  an  oscillating  airfoil 
using"  unsteady  methods  (16)/  (17) ,  have  demonstrated  good  quali¬ 
tative  agreement  with  experimental  observations.  However/  one 
feature  at  least  needs  further  examination  in  regard  to  improved 
modeling:  it  is  reported  (18)/  (19)  that  when  incidence  is 
increasing  beyond  the  static  stall  angle/  the  location  of  zero 
skin  friction  in  the  turbulent  boundary  layer  and  the  catas¬ 
trophic  separation  can  occur  at  different  stations/  Figure  1. 
Apparently,  a  long  thin  tongue  of  this  reversed  flow  precedes  the 
main  separation  zone.  This  is  not  observed  under  quasi-steady 
conditions. 

Crimi  and  Reeves  (20)  combined  a  potential  flow  method  with 
an  unsteady  boundary  layer  analysis  in  a  viscous/inviscid  inter¬ 
action  approach.  The  potential  flow  model  was  based  on  chordline 
singularities  and  so  excluded  modeling  of  the  thick  wake.  Alsof 
detail  of  the  stagnation  point  location  in  relation  to  the  curved 
leading  edge  was  missing.  Emphasis  was  placed  on  the  details  of 
leading-edge  bubble  bursting  and  application  to  trailing-edge 
separation  was  not  attempted.  Shamroth  and  Kreskovsky  (21)  de¬ 
veloped  a  similar  technique  but  with  improved  treatment  of  the 
separated  flow  region,  transition  phenomena  and  potential  flow 
region.  However,  the  procedure  failed  to  predict  the  flow  field 
about  the  stalled  airfoil.  They  concluded  that  the  effect  on  the 
outer  inviscid  solution  due  to  the  finite  wake  displacement  must 
be  modelled. 

In  spite  of  the  shortcomings  of  the  above  approaches,  the 
general  technique  of  matching  various  viscous  and  inviscid  re¬ 
gions  remains  an  attractive  alternative  to  the  full  Navier-Stokes 
treatment.  Although,  in  principle,  the  latter  can  overcome  limi¬ 
tations  of  the  potential  flow/boundary  layer  iterative  approach. 
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Figure  1.  The  Model  of  Sears  and  Telionis  (19)  for 
Upstream-Moving  Separation  of  an  Unsteady 
Boundary  Layer , 


such  treatment  is  limited  at  this  time  to  laminar  flows  at 
Reynolds  numbers  much  lower  than  realistic  for  most  practical 
applications  (2).  Progress  towards  higher  Reynolds  numbers  is 
being  made,  but  applications  to  general  problems  is  still  a  long 
way  away  (22). 

The  present  method  goes  beyond  the  capabilities  of  the 
earlier  theoretical  approaches  in  that  both  trailing-edge  and 
leading-edge  stall  with  vortex  passage  can  be  included  in  prin¬ 
ciple.  The  method#  developed  for  the  three-dimensional  case#  is 
applicable  to  arbitrary  configurations  and  to  general  motions; 
i.e.#  not  just  pitch  oscillation.  In  addition#  because  the  basis 
of  the  method  is  a  surface  singularity  panel  code#  a  more  reli¬ 
able  and  direct  coupling  between  the  inviscid  and  viscous  analy¬ 
sis  is  assured.  Moreover#  modeling  of  the  separated  zone  in  the 
trailing-edge  region  and  more  detailed  treatment  of  the 
vortex/surface  interaction  should  make  the  approach  more  viable 
for  applications  to  dynamic  separation  problems. 


2.1  General 


Although  remarkable  advances  are  being  made  in  flow  field 
calculations  using  finite-difference  and  also  finite-element 
methods#  the  surface  integral  approach  using  panel  methods 
coupled  with  special  routines  for  nonlinear  effects  still  offers 
distinct  advantages  for  many  real  flow  problems.  In  particular, 

ranel  methods  offer  greater  versatility  for  practical  application 
o  complicated  configurations  and  are  considerably  more  efficient 
ir-  terms  of  computing  effort.  However#  the  concept  of  zonal 
modeling — in  which  a  local  Navier-Stokes  analysis  is  coupled  with 
a  panel  method — should  not  be  overlooked.  Ultimately#  such  a 
coupling  should  lead  to  an  improved  modeling  of  vortices  (e.g., 
vortex  cores#  vortex  dissipation  and  breakdown)#  thick  viscous 
regions#  local  separation  bubbles#  and  shock  wave/boundary  in¬ 
teractions. 

Over  the  past  decade#  panel  methods  have  seen  a  trend  toward 
higher-order  formulations  (23),  (25)  and  (26).  At  the  outset  it 
was  argued  that  compared  with  the  earlier  low-order  methods  the 
more  continuous  representation  of  the  surface  singularity  distri¬ 
bution  should  allow  a  reduction  in  panel  density  for  a  given 
solution  accuracy  and#  hence#  should  lead  to  lower  computing 
costs.  No  such  benefits  have  appeared  so  far  for  the  general 
three-dimensional  case.  In  fact#  preliminary  investigations  (27) 
have  indicated  that  the  prediction  accuracy  for  problems  with 
complicated  interactions#  such  as  vortex/surface  or  high  curva¬ 
ture  situations#  depends  more  on  the  density  of  control  points 
where  the  boundary  conditions  are  enforced;  the  order  of  the 
singularity  distribution  plays  a  minor  role.  In  the  meantime# 
further  developments  of  piecewise  constant  singularity  panel 
methods,  e.g.,  Morino  (28)  and  AMI's  Program  VSAERO  (6),  are 
giving  comparable  accuracy  to  the  higher  order  methods  at  much 
lower  computing  costs. 

The  low  computing  cost  of  Program  VSAERO  makes  it  practical 
to  apply  the  method  to  nonlinear  problems  requiring  iterative 
solutions#  e.g.,  wake  relaxation  for  high-lift  configurations, 
multiple-component  problems  and  rotor  cases;  and  viscous/inviscid 
calculations  with  coupled  boundary  layer  analyses#  including  the 
case  with  extensive  separation  (5),  and  also  time-stepping  calcu¬ 
lations  for  three-dimensional  unsteady  problems  (4)  that  are 
beyond  the  scope  of  a  harmonic  analysis.  This  method,  therefore, 
offers  an  attractive  basis  for  a  practical  tool  aimed  at  pre¬ 
dicting  the  aerodynamic  characteristics  of  dynamic  stall  prob¬ 
lems.  At  the  outset,  this  code  was  being  developed  in  two  dif¬ 
ferent  directions;  viz,  one  was  for  extensive  separation  modeling 
under  the  "steady"  conditions,  while  the  other  was  for  unsteady 
time-stepping  calculations.  These  two  capabilities,  described  in 
the  section  below#  have  now  been  brought  together  in  one  code. 


Under  essentially  steady  conditions,  the  pressure  distribu¬ 
tion  in  a  trailing-edge  separation  region  is  usually  charac¬ 
terized  by  a  constant  pressure  region  extending  back  to  the 
trailing  edge  followed  by  a  short  recompression  region  le.g., 
(29)).  A  simplification  of  this  characteristic  is  modelled  in 
the  two-dimensional  CLMAX  program  (30),  (31)  using  a  pair  of 
constant  strength  vortex  sheets  to  enclose  the  separated  region. 
Figure  2.  The  length  of  the  sheets  required  a  semi-empirical 
approach  and  the  condition  that  the  sheets  be  force-free  is 
satisfied  in  an  interactive  cycle  in  which  segments  of  the  sheet 
are  aligned  with  calculated  local  flow  directions.  The  method 
combines  boundary  layer  and  potential  flow  codes  in  an  outer 
inviscid/viscid  iteration  cycle.  The  potential  flow  pressure 
distribution — which  includes  the  influence  of  the  free  vortex 
sheets — is  passed  to  the  boundary  layer  routine  which  then  sup¬ 
plies  the  separation  points  and  the  boundary  layer  displacement 
thickness  distribution  for  the  next  iteration.  The  boundary 
layer  displacement  effect  in  the  attached  flow  region  is  modelled 
by  transpiration  (i.e.,  source  distribution)  rather  than  by  a 
displacement  surface.  The  main  advantage  of  the  transpiration 
approach  is  that  the  matrix  of  influence  coefficients  in  the 
panel  method  remains  essentially  the  same  from  one  iteration  to 
the  next;  only  the  wake  condition  changes. 

The  thin  vortex  sheet  model  of  the  upper  separated  shear 
layer  was  demonstrated  by  Young  and  Hoad  (32)  to  be  a  reasonable 
representation  of  the  flow  as  far  back  as  the  trailing  edge.  For 
example,  a  comparison  from  (32)  of  a  laser-velocimeter  flow 
survey,  and  CL MAX  program  calculation  1b  shown  here  in  Figure  3. 
Downstream  of  the  trailing  edge  the  vortex  sheet  model  becomes 
less  representative  of  the  flow;  however,  a  later  evaluation  of  a 
graded  vorticity  model  over  the  recompression  zone  showed  little 
effect  on  the  airfoil  solution.  More  detailed  modeling  of  the 
recompression  zone  (such  as,  for  example,  the  approaches  used  by 
Gross  (33)  or  Zumwalt  (34)  would  be  desirable  in  cases  where  the 
wake  interacts  closely  with  a  downstream  component. 

A  particular  feature  of  the  vortex  sheet  model  enclosing  the 
region  of  low  energy  is  that  pressure  can  be  calculated  directly 
in  the  separated  zone  (30).  This  is  an  additional  advantage  over 
the  displacement  surface  approach  of  Henderson  (35)  and  over  the 
source  outflow  model  of  Jacob  (36).  The  CLMAX  method  generally 
gives  very  close  agreement  with  experimental  pressure  distribu¬ 
tions  (30),  (31).  An  initial  extension  of  this  method  to  the 
unsteady  case  is  reported  in  (31)  where  quasi-steady  solutions 
coupled  with  a  phase  shift  model  were  used.  Extension  of  the 
model  to  the  three-dimensional  case  is  reported  in  (37)  for  a 
stripwise  model  and  in  (5)  for  a  more  general  treatment.  The 
separation  model  has  also  been  successfully  installed  in  a  tran¬ 
sonic  finite-difference  code  (38). 
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Figure  3.  Location  of  Experimental  Free  Shear  Layer  and 
CLMAX  Calculated  Vortex  Sheet  Centerlines  for 
Low  Mach  Number  Case  (from  (32)). 


The  sane  basic  nodel  is  also  applicable  for  the  unsteady 
case;  here,  however#  the  assumption  of  constant  vorticity  is  no 
longer  valid.  In  fact#  a  dynamic  wake  model  is  essential  and 
will  be  discussed  below  after  the  description  of  the  unsteady 
formulation. 


2.3  Unsteady  Method 

2.3.1  formulation 


Consider  the  whole  of  space  divided  into  two  regions  by  the 
surface  of  the  configuration  and  assume  the  existence  of 
Laplacian  velocity  potential  distributions  in  the  two  regions; 
i.e.#  4>  in  the  flow  field  and  * i  in  the  blade  interior.  If  we 
now  apply  Green's  third  identity  to  the  two  regions#  then  the 
total  potential  at  a  point#  P#  on  the  inside  surface  of  the 
boundary  can  be  written: 


//  (♦-•,)»•  <?) « -  H*  ■ 

S-P 


//  (‘u  -  *l)“  •  *(f)  dU 


■p  ^dS  +  4*4 


(1) 


Here,  r  is  the  length  of  the  vector  from  the  surface  element  to 
the  point#  P,  and  S-P  signifies  that  the  point#  P#  is  excluded 
from  the  surface  integration.  Equation  (1)  includes  the  contri¬ 
bution  from  the  wake  surface#  W. 


The  Dirichlet  boundary  condition  is  now  applied  in  the 
interior  region  to  render  a  unique  distribution.  In  principle# 
any  potential  flow  can  be  applied.  However#  the  flow,  *i  ■  ^oo# 
implied  by  Morino  (28)  and  used  by  Johnson  and  Rubbert  (25),  has 


proven  to  be 
becomes 

very  reliable  in  practice. 

With  this  flow#  Eq 
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where  <p ,  the  perturbation  potential  in  the  flow  field,  has  been 

substituted  for  $  -  4>  . 

00 

The  first  two  terms  in  Eq.  (2)  give  the  perturbation  poten¬ 
tial  due  to  a  distribution  of  normal  doublets  of  strength,  <t>  ,  on 
the  configuration  surface.  Similarly,  the  third  term  represents 
a  doublet  distribution  of  strength,  fyj  -  <t>L  ,  on  the  wake  and  the 
fourth  term  represents  a  source  distribution  of  strength,  n  •  , 
on  the  configuration  surface. 

Equation  (2)  is  basically  the  same  as  the  formulation  given 
by  Morino  (28)  who  used  a  direct  application  of  Green's  theorem 
in  the  flow  field.  The  present  approach  to  the  problem  is  a 
special  case  of  multi-domain  formulation  which  has  led  to  the 
more  general  three-dimensional  method  in  which  large  regions  of 
separated  flow  are  modelled  in  a  similar  way  to  that  in  tne  CLMAX 
program  (5),  (30)  and  (31). 

The  source  term  in  Eq.  (2)  can  be  evaluated  directly  from 
the  condition  of  no  flow  penetration  at  the  surface.  The  flow 
velocity  relative  to  the  body-fixed  frame  is  at  any  instant  of 
time 


V*v  +  V00-fth~  R  (3) 

where  the  perturbation  velocity  is  v  h  is  the  axis  of 

rotation,  v  «,  and  &  are  the  instantaneous  onset  flow  and  angular 
velocity,  respectively,  and  R  is  the  relative  position  vector 
between  a  point  on  the  rotation  axis  and  a  point  on  the  surface. 

For  zero  penetration,  V  •  n  ■  0.  Hence, 


n  •  V$  ■  n  •  Vw  -  8n  •  b  *  R, 


and  Eq.  (2)  becomes 
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This  is  the  basic  equation  of  the  method.  It  is  solved  for 
the  unknown  surface  perturbation  potential,  ♦  ,  or  surface  doublet 
distribution,  u,  at  a  number  of  time  steps  as  the  configuration 
proceeds  through  the  motion.  The  wake  surface  is  transported  at 
the  end  of  each  step  using  calculated  velocities  of  points  on  the 
wake  surface.  The  doublet  distribution,  4*u  -  ♦l,  on  each  wake 
surface  is  known  from  solutions  at  earlier  time  steps.  The 
unsteady  Kutta  condition 


In.  +  v  -^  =  o 

3t  V  3s 


(5) 


is  satisfied  at  points  along  each  wake  separation  line  at  each 
time  step. 

At  each  time  step  the  flow  solution  is  determined  with 
reference  to  the  body-fixed  frame.  The  incompressible  pressure 
coefficient  is,  therefore,  given  by 


cP  ■  (v  - v!  +  2  UK 


(6) 


where  Vs  «  fih  is  the  instantaneous  velocity  of  a  point 
on  the  surface  relative  to  a  stationary  reference  frame,  and  V  is 
given  by  Eq.  (3). 

2.3.2  Mumerlcal  Procedure 

The  general  arrangement  of  the  configuration  is  Bhown  in 
Figure  4.  The  x,y,z  coordinate  system  with  unit  vectors,  i,  j, 
k,  is  fixed  relative  to  the  configuration.  For  symmetrical 
applications,  the  x-z  plane  is  regarded  bb  the  plane  of  symmetry. 

A  numerical  procedure  has  been  assembled  in  a  time-stepping 
mode  to  obtain  the  unsteady  pressure  distribution  and  forces  and 
moments.  The  surface  of  the  wing  is  represented  by  planar  quad¬ 
rilateral  panels  over  each  of  which  the  doublet  and  source  dis¬ 
tributions  are  assumed  constant.  With  this  assumption,  the  sur¬ 
face  integrals  in  Eq.  (4)  can  be  performed  in  the  closed  form  for 
each  panel. 

Equation  (4)  is  then  satisfied  simultaneously  at  a  point  in 
the  center  of  each  panel.  If  there  are  N  panels  representing  the 
configuration  surface,  Eq.  (4)  becomes: 


BODY-FIXED  COORDINATE  SYSTEM 


where  pK  is  the  unknown  doublet  value  on  Panel  K.  (Note:  pt 
4>k/4tt.)  *  1 


If 

1 J  =  ^  I  ^WK  CJK  f  +  'V.  *  aH 


where  N^,  the  number  of  panels  in  the  wake,  varies  with  time  and 
Q  and  V  take  their  instantaneous  values  at  each  time  step. 


*  \  bjkI/4* 


is  the  source  distribution  due  to  rotation  about  the  axis,  h,  and 


V  '  £  i  “K  BJK  j /4” 


are  the  components  of  a  three-part  source  distribution  due  to  the 
relative  translation  of  the  configuration  and  the  onset  flow. 
(Note:  in  a  symmetrical  case  the  y-component  is  zero.) 

The  quantities,  Bjk  and  Cjk,  are  the  velocity  potential 
influence  coefficients  for  the  constant  source  and  doublet  dis¬ 
tributions,  respectively,  on  panel  K  acting  on  the  control  point 
on  panel  J.  These  include  contributions  from  the  image  panel  in 
the  symmetrical  case.  Expressions  for  these  influence  coeffi¬ 
cients  have  been  given  by  Morino  in  (28)  based  on  hyperbolic 
paraboloidal  panels.  Slightly  different  expressions  are  instal¬ 
led  in  the  VSAERO  code  based  on  planar  panels. 

Equation  (7)  is  solved  by  a  direct  method  for  N  i  320  and  by 
an  iterative  method  for  N  >  320. 

The  surface  pressure  distribution  is  calculated  using  Eq. 
(6).  The  surface  gradient  of  p  is  evaluated  on  each  panel  by 
differentiating  a  two-way  parabolic  fit  through  the  doublet 
values  on  the  panel  and  its  four  immediate  neighbors.  At  the 
separation  lines  a  simple  differencing  is  applied  for  the 
gradients  approaching  the  separation  line. 

The  gradient  of  <P  with  respect  to  time  is  evaluated  by 
central  differencing  over  two  time  steps;  i.e.. 


For  harmonic  motions  the  real  and  imaginary  pressure  are 
obtained  by  Fourier  analysis  for  the  first  harmonic  based  on 
solutions  over  a  complete  cycle.  The  calculations  start  with 
incidence  a0  and  a  regular  (i.e.,  steady)  wake.  Two  iterations 
are  performed  to  render  the  wake  force  free.  An  oscillatory 
doublet  component  based  on  a  linearized  solution  is  then  super¬ 
imposed  along  each  wake  line  before  starting  the  time-step  model. 
Time-step  calculations  proceed  over  a  half  cycle  before  applying 
the  Fourier  analysis. 

At  each  time  step  a  new  panel  is  formed  at  the  head  of  each 
column  of  wake  panels  and  all  the  existing  wake  panel  corner 
points  are  convected  downstream  at  the  local  velocity.  Each  wake 
panel  keeps  the  doublet  value  it  received  at  the  time  it  was 
formed.  This  doublet  value  is  based  on  the  conditions  at  the 
separation  line  and  satisfies  Eq.  (5).  It  is  assumed  that  the 
shedding  occurs  at  constant  vorticity  over  the  time  interval#  t. 


In  this  way  the  doublet  strength#  u 


t+At 


W 


#  on  the  new  wake  panel 


is  related  to  the  strength#  uw  #  of  the  previous  wake  panel  at 
the  separation  line  by 


t+At  _  -  t  t 


where  uT  is  the  resultant  doublet  value  at  the  separation  line. 


2.4  Unsteady  Separated  Flow 

The  combined  code  for  separated  flow  modeling  in  the  un¬ 
steady  case  requires  a  more  sophisticated  treatment  of  the  free- 
shear  layer  model  than  was  used  for  the  steady  case.  Velocities 
are  still  calculated  at  a  set  of  points  along  each  free  sheet# 
but  in  the  unsteady  case  we  now  transport  these  points  (and  their 
associated  doublet  value)  along  the  calculated  velocity  vectors 
for  a  small  time  interval#  At.  In  this  way  as  time  progresses#  a 
dynamic  wake  model  is  generated.  At  each  step  a  new  piece  of 
free  sheet  is  shed  from  the  calculated  separation  point;  the 
strength  and  size  of  this  new  segment  is  determined  by  the  local 
upstream  velocity  condition.  The  location  of  the  separation — 
calculated  using  an  unsteady  boundary  layer  code#  see  the  next 
section — can  now  move  with  time. 


It  is  convenient  to  regard  the  local  vorticity  (i.e.#  doub¬ 
let  gradient  on  the  free  vortex  sheets  in  two  components;  a 
streamwise  component  and  a  cross-flow  component.  The  streamwise 
component  is  already  force  free  and  is  related  to  the  spanwise 
rate  of  shedding  of  circulation  from  the  wing.  The  cross-flow 
vorticity  component  is  associated  with  direct  dumping  of  bound 


circulation  from  the  separation  line  and  must  be  transported  with 
the  local  flow  velocity  in  order  to  be  force  free.  This  cross¬ 
flow  vorticity  component— which  was  assumed  constant  with  stream- 
wise  distance  in  the  steady  case — now  varies  along  each  stream¬ 
line  on  the  free  sheets  for  two  reasons:  first,  the  vorticity 
value  being  convected  onto  the  free  sheet  at  each  separation 
point  is  varying  in  time  because  of  varying  onset  flow  conditions 
and  because  of  the  changing  separation  locations;  secondly, 
stretching  by  the  entire  configuration  of  solid  surfaces  and  free 
wake  sheets.  This  stretching  of  the  doublet  distribution  carried 
by  the  free  sheets  yields  varying  vorticity  values  when  the 
doublet  gradient  is  evaluated.  In  this  way  the  free  sheets  can 
become  highly  distorted  and  centers  of  vortex  roll-up  may  form. 
Special  treatment  of  the  sheets  is  therefore  essential  if  numeri¬ 
cal  stability  is  to  be  maintained.  Two  routines  are  being 
evaluated  in  this  work  but  they  have  not  been  fully  implemented 
at  this  time. 

(i)  Vortex  Amalgamation 

To  cater  for  vortex  roll-up  in  a  reasonable  manner  it  is 
essential  to  include  a  vortex  core  model  in  which  integrated 
vorticity  is  accumulated  rather  than  to  follow  a  detailed  calcu¬ 
lation  of  multiple  turns  of  a  vortex  spiral.  An  amalgamation 
scheme  similar  to  that  of  Moore  (39)  is  being  used.  When  the 
angle  between  neighboring  segments  representing  a  sheet  exceeds  a 
specified  angle,  the  segment  end  points  are  merged  to  a  new 
location  at  the  centroid  of  their  combined  circulation.  A  number 
of  such  cores  are  allowed  in  the  new  routine  to  deal  with  com¬ 
plex  motions.  A  viscous  core  expression  can  be  applied  to  each 
vortex  core  when  computing  the  field  velocities. 

(ii)  Redistribution 

Having  performed  the  vortex-core  amalgamation  calculation 
along  each  free  sheet  the  points  defining  the  intermediate  free 
sheets  are  redistributed  with  equal  spacing  in  a  manner  similar 
to  Fink  and  Soh  (40)  and  Sarpkaya  and  Schoaff  (41).  Portions  of 
the  sheet  between  amalgamated  cores  are  treated  independently 
here.  Figure  5.  This  treatment,  which  is  applied  to  both  the 
sheet  geometry  and  its  doublet  distribution,  uses  a  biquadratic 
interpolation  scheme  based  on  surface  distance  along  the  sheet. 

This  routine  should  help  stabilize  the  numerical  calcula¬ 
tions,  especially  in  the  initial  part  of  each  sheet  when  the 
separation  location  is  varying  with  time.  In  the  three-dimen¬ 
sional  case  the  redistribution  scheme  is  being  arranged  along  the 
calculated  mean  streamlines  in  the  wake  sheets.  This  causes  some 
difficulty  if  amalgamation  is  not  proceeding  uniformly  along  all 
lines  on  a  sheet. 


2.5 


As  the  routines  are  being  developed,  preliminary  calcula¬ 
tions  are  being  performed  to  check  the  basic  operation  of  the 
code.  Figure  6  shows  the  growth  of  indicial  lift  and  circulation 
for  a  NACA  0012  impulsively  started  from  rest  at  an  angle  of 
attack  of  .1  rad.  The  curves  are  compared  with  Wagner's  function 
for  indicial  lift  and  R.T.  Jones  indicial  circulation  for  a  flat 
plate.  The  initial  calculations,  which  used  31  panels  around  the 
section,  are  in  good  agreement  and  indicate  a  slightly  higher 
trend  which  is  consistent  with  a  higher  steady  state  circulation 
for  the  thickness  case. 

Some  recent  refinements  developed  in  the  two-dimensional 
pilot  code  have  significantly  reduced  the  computing  requirement 
of  these  time-stepping  calculations.  For  example.  Figure  7  shows 
the  effect  on  indicial  lift  of  varying  the  number  of  time  steps 
in  the  Wagner  problem  and  demonstrates  a  rapid  convergence. 

The  procedure  has  been  tested  also  for  the  harmonic  oscilla¬ 
tion  case.  Earlier  calculations  required  80  time  steps  per  cycle 
for  a  NACA  0012  oscillating  in  pitch  about  the  quarter  chord. 
These  compared  favorably  with  the  Theodorsen  flat  plate  function 
over  a  range  of  reduced  frequency.  Figure  8.  The  new  calcula¬ 
tions  are  also  in  good  agreement  but  were  performed  with  only  16 
time  steps  per  cycle.  Figure  9  shows  the  computed  results.  Cl 
versus  time,  using  only  4  time  steps  per  cycle.  This  is  in 
remarkably  good  agreement  with  the  16  and  also  32  time-step/cycle 
solutions,  demonstrating  an  extremely  good  convergence  charac¬ 
teristic. 

Time-stepping  calculations  have  also  been  performed  for 
cases  with  prescribed  extensive  separations.  The  purpose  of 
these  calculations  was  to  check  the  basic  unsteady  circulation 
shedding  model  in  the  potential  flow  code.  For  the  first  set  of 
tests,  the  wake  panels  were  simply  transported  at  the  onset  flow 
velocity  after  the  initial  growth  as  determined  from  the  surface 
conditions  at  separation.  Several  triangular  shapes  were  con¬ 
sidered,  each  starting  impulsively  from  rest  and  proceeding  for¬ 
ward  over  10  time  steps  for  a  total  time  of  t  *=  x  U^/h  ■  3.0, 
where  h  is  the  triangle  base  height.  Separation  was  prescribed 
at  the  corners.  Figure  10(a)  shows  the  computed  history  of  the 
drag  coefficient  from  pressure  integration  for  a  60  degree  tri¬ 
angle  with  blunt  face  forward.  A  total  of  40  panels  was  used  to 
represent  the  triangle  surface.  The  calculation  was  repeated  in 
the  presence  of  wind  tunnel  walls  (also  panelled)  with  a  10% 
blockage  ratio.  The  indicated  blockage  correction  is  somewhat 
lower  than  that  given  by  standard  techniques.  Figure  10(b) 
compares  the  computed  pressure  distributions  for  this  triangle  in 
and  out  of  the  tunnel.  This  "base"  pressure  has  only  a  small 
variation  and  is  quite  close  to  experimental  measurements. 
Figure  11  shows  a  summary  of  computed  drag  coefficient  versus 
triangle  semi-apex  angle.  The  calculated  values  are  slightly 
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8.  Comparison  of  Real  and  Imaginary  Lifts  as  a 

Function  of  Reduced  Frequency,  a  =  1°  sin  kt. 


STARTED  IMPULSIVELY  FROM  REST:  c  -  U|/h=3.0 


high  in  relation  to  the  experimental  data  collected  from  several 
sources  by  Hoerner  in  Fluid  Dynamics. 

One  further  case  was  run  for  the  60  degree,  apex-forward 
triangle  in  free  air  with  the  full  wake  velocity  calculation 
routine  turned  on  but  without  the  amalgamation  and  redistribution 
schemes  at  this  stage.  The  calculated  CD  for  this  case  falls 
below  the  experimental  value.  Figure  11.  A  series  of  computed 
wake  shapes  is  shown  in  Figure  12.  These  are  samples  from  a 
total  of  40  time-step  calculations.  The  total  computing  time  for 
this  case  was  195  seconds  on  a  PRIME  550  minicomputer — this  is 
equivalent  to  less  than  2  seconds  of  CRAY  time.  The  solution 
should  benefit  from  the  numerical  damping  provided  by  the  amalga¬ 
mation  and  redistribution  schemes  described  earlier. 

Finally,  a  test  calculation  was  performed  for  a  NACA  0012 
section  in  a  state  of  pitch  from  10  degrees  to  30  degrees  with 
ac/240>  .175.  The  calculation  used  30  panels  and  10  time  steps. 
Separation  points  were  prescribed  and  the  motion  was  started 
impulsively  from  rest. 

Figure  13(a)  shows  a  sample  of  the  computed  wake  shapes  and 
demonstrates  a  reasonable  numerical  behavior.  Sample  pressure 
distributions  are  shown  in  Figure  13(b).  The  passage  of  the 
leading-edge  vortex  is  clearly  shown.  This  is  associated  with  a 
local  region  of  reversed  flow.  These  are  preliminary  test  calcu¬ 
lations  aimed  at  exploring  the  numerical  behavior  of  the  calcula¬ 
tion  procedure  and  potential  flow  model.  Future  cases  will 
include  the  coupled  boundary  layer  calculation  for  predicting  the 
separation  point.  At  that  time  the  calculated  results  will  be 
compared  with  experimental  data. 
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Figure  13(a).  NACA  0012  STARTING  IMPULSIVELY  FROM  REST 

PITCHING  FROM  10°  TO  30°  :  «r/o n  «, 


NACA  0012  STARTING  IMPULSIVELY  FROM  REST 
PITCHING  FROM  10°  TO  30°:  d  C/2U.=0.175 


Figure  13(b)  (ii). 


Curie's  original  method  has  been  modified  to  calculate  the 
unsteady  boundary  layer  development.  This  is  achieved  by  solving 
the  unsteady  momentum  integral  equation  using  a  Runge-Kutta 
method.  The  turbulent  boundary  layer  method  is  based  on  the 
unsteady  momentum  integral  equation  as  in  the  laminar  boundary 
layer  method.  Cousteix's  entrainment  relationship  (43)  and 
Lyrio/Ferziger's  skin  friction  relationship  (44)  are  used  for 
closure.  The  details  of  the  methods  are  described  in  the  fol¬ 
lowing  sections. 


ue2  at 


This  is  a  first-order  differential  equation  for  V  and  can  be 
solved  for  9  if  additional  relationships  for  6 ,  Cf  and  <$*  are 
given.  In  Curie's  steady  boundary  layer  procedure, rthese  rela¬ 
tionships  are  expressed  as 


Cf  =  2  •  £/R0 


L  =  2 [£  -  K(H  +  2)  ] 


where 


!»,£  ■  functions  of  u  and  K 


K  -  9  /v (dUe/dx) 


v  -  K2UeU^/U4 


Calculation  begins  at  the  stagnation  point  and  K  takes  the 
starting  value,  K0  *  0.0855.  The  inital  momentum  thickness,  0O  > 


(0.0855V/ (dU^/d  ) 


For  unsteady  turbulent  boundary  layers,  the  momentum  inte 
gral  equation  and  the  entrainment  equation  are  given  by: 


U  3xl~e 


where 


C  _  ii  Ze 

'"E  3x  U  * 

e 


Equations  (14}  and  (15)  have  five  unknowns;  i.e.,  C &  Cf,  6,6, 
and  6  *,  and  the  system  needs  additional  relationships  for 
closure. 

For  the  skin  friction,  a  new  correlation  of  Lyrio  et  al. 
(44)  is  used: 


~  =  .051(1  -  2A 


1.732  0.268 

(A/Re*)  sgn (1  -  2A) 


where 


=  1.5A  +  0.279Vt  +  0.321  V*/A 


A  =  ii 
A  6 


Ji s  i  sgn(Tw) 


and  k  *  0.41  is  the  von  Karroan  constant. 


The  similarity  solutions  have  shown  that  the  entrainment 
coefficient/  CE,  can  be  expressed  as 

C  =  c  _  JL  ££ 
e  ces  ue  at 

(21) 

where  Ct,  is  the  entrainment  coefficient  for  the  steady  case  as 
ES 

given  by: 


C 


Es 


(0.074G  -  1.0957/G) 


(22) 


Substitution  of  Eq.  (21)  into  Eqs.  (14)  and  (15)  and  the 
introduction  of  H*  *  (6  -  6*)/ 8  give 


_L  Ml  4.  li  =  T» 

U  3  t  3x  1 
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(23) 
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(24) 


with 


Cf  8  (H  +  2)  8ue 

B1  "  2  Ue  3x 
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Equations  (23)  and  (24)  constitute  a  system  of  first-order 
hyperbolic  partial  differential  equations  for  6*  and  9.  The 
characteristic  equation  of  the  second  deqree  for  the  parameter,  * 

dx 
e 


X2(H*  -  HH*')  +  X [1  -  H*  +  H*' (1  +  H) )  -  H*'  -  0 


(25) 


Equations  (23)  and  (24)  can  be  solved  in  various  ways. 
Initial  conditions  at  t  *  5o  and  boundary  conditions  at  the 
stagnation  point  (t  >  0)  are  sufficient  to  determine  a  solution 
in  the  region  where  the  flow  is  attached;  i.e.,  *2  >  0*  In  the 
present  paper,  the  time  derivatives  are  treated  as  forcing  terms 
and  the  integration  is  performed  in  the  x-direction  using  a 
Runge-Kutta  method. 

The  boundary  layer  procedure  has  been  tested  against 
experiments  and  the  calculations  of  other  investigators.  The 
results  are  shown  in  Figures  14  and  15.  Figure  14(a),  (b),  and 
(c)  shows  the  mean  quantities  (momentun  thickness,  skin  friction 
and  shape  factor)  for  the  experiment  conducted  by  Cousteix  (45) 
on  a  flat  plate.  The  main  free  stream  velocity  is  22  m/sec  and 
the  motion  is  harmonic  with  respect  to  time  with  the  frequency  of 
38  hz.  The  present  calculation  (solid  line)  agrees  very  well 
with  the  calculation  by  Cousteix  and  the  correlation  with 
experiment  is  also  very  good. 

Figure  15  shows  the  comparison  with  the  calculation  of  Nash 
et  al.  (46)  for  a  monotonically  time-varying  flow  on  a  flat 
plate.  The  present  calculation  predicts  the  separation  at  the 
end  of  the  plate  when  wt  =  0.682  as  in  the  Nash  et  al. 
calculation,  the  overall  results  are  in  good  agreement  with 
their  calculation. 

These  test  calculations  have  shown  that  the  current  boundary 
layer  procedure  is  capable  of  predicting  the  unsteady  boundary 
layer  development  and  it  should  be  adequate  for  analyzing  dynamic 
stall  problems. 


The  complete  procedure  coupling  the  unsteady  time-stepping 
potential  flow  panel  method*  the  extensive  separation  wake  model 
and  the  unsteady  boundary  layer  code  has  been  assembled  in  the 

Silot  code  for  a  system  checkout  prior  to  forming  the  three- 
imensional  version.  The  flow  diagram  for  the  procedure  is  shown 
in  Figure  16.  At  this  time  the  unsteady  boundary  layer  code  is 
called  at  each  time  step  and  is  fed  by  unsteady  derivatives  from 
the  potential  flow  calculation. 

An  experimental  data  case  from  (3)  was  run  and  the  computed 
lift  variation  with  a  compared  with  the  measured  data  in  Figure 
17.  The  airfoil  is  a  NACA  0012  and  is  oscillating  in  pitch  about 
the  quarter-chord  line  with  a  »  8.1°  +  4.9°  sin  (0.2t);  i.e.* 
below  the  dynamic  stall  onset.  Reynolds  number  was  4  x  10°. 
This  reduced  frequency  condition  is  very  close  to  the  changeover 
from  a  lead  to  a  lag  situation  and  so  there  is  only  a  small 
difference  between  the  upswing  and  downswing  curves. 

A  preliminary  calculation  was  performed  for  a  NACA  0012 
section  oscillating  in  pitch  about  the  quarter-chord  with  a  «  10° 
sin  (2t).  Figure  18  shows  the  predicted  history  of  the  separa¬ 
tion  location  superimposed  on  the  a  history.  The  most  forward 
separation  reached  .2  x/c  with  a  phase  lag  of  about  17*. 

The  above  calculation  was  performed  with  just  the  trailing- 
edge  wake.  Calculations  are  now  proceeding  with  linear  pitch-up 
with  a  second  (upper  surface)  separated  wake.  Experimental 
measurements  of  airfoils  undergoing  constant  rate  pitch-up 
motions  have  been  taken  at  the  Frank  J.  Seiler  Research  Labora¬ 
tory  (47).  Three  cases  are  considered  briefly  here  for  correla¬ 
tion  purposes.  In  these  cases  the  airfoil  is  pitched  up  from  a  ■ 
0  to  approximately  1  radian  at  a  constant  pitch  rate  and  then 
held  at  constant  angle  of  attack.  Three  j>itch  rates  are  con¬ 
sidered  with  normalized  pitch  rates*  k  (*  a  C/2U,*,)  of  .047*  .089 
and  .133.  Figures  19  (a),  (b),  and  (c)*  respectively*  show  the 
comparison  between  calculated  and  measured  CL  ~  <*  and  CD  ~  <* 
characteristics  for  these  cases.  For  the  low  pitch  rate*  k  * 
0.047*  the  comparison  is  very  good  up  to  about  30®  but  then 
deteriorates.  The  calculated  Cl  remains  fairly  constant  with  « 
until  the  pitch  rate  drops  to  zero  while  the  experimental  curve 
falls  markedly.  The  calculated  rise  in  drag  with  a  has  a  steady 
rate  in  the  30®  to  60°  range  while  the  experimental  measurements 
include  a  substantial  increment  above  this  rate  peaking  at  about 
a  *  40o.  There  is  a  good  agreement  between  calculated  and 
measured  lift  and  drag  values  for  the  final  "steady  state"  condi¬ 
tions  at  a  -  600. 

One  possible  reason  for  the  departure  of  the  calculated  lift 
and  drag  in  the  latter  part  of  the  pitch-up  phase  is  the  modeling 
of  the  leading-edge  vortex  roll-up.  The  amalgamation  and  redis¬ 
tribution  schemes  that  were  installed  to  stabilize  the  dynamic 
wake  calculations  are  not  performing  in  a  consistent  manner  at 


Figure  17 .  Comparison  of  Calculated  and  Measured  Lift 
on  a  NACA  0012  Airfoil  Oscillating  in  Pitch 
about  the  Quarter  Chord. 
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Figure  18.  Calculated  History  of  Separation  Location  on  a  NACA  0012  Airfoil 
Oscillating  about  the  Quarter-Chord. 


O  CALCULATION 
— —  EXPERIMENT  (47) 


(a)  k  =  0.047,  =  60°. 

Figure  19.  Comparison  of  Calculated  and  Measured  Lift  and 
Pressure  Drag  on  a  NACA  0012  Section  During 
Pitch-up  Motion  about  x/c  =  .317. 
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c)  k  =  0.13,  cij^  =  55°. 
Figure  19.  Concluded. 
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this  time  and  so  further  refinements  are  planned.  In  the  present 
calculation  the  vortex  core  on  the  upper  sheet  did  not  "condense" 
early  enough  in  the  amalgamation  procedure;  consequently,  the 
vortex  formed  just  downstream  of  the  trailing  edge  and  did  not 
closely  interact  with  the  airfoil  surface  during  the  pitch-up 
phase.  This  tendency  was  still  present  at  the  higher  pitch-up 
rates,  k  *  .089  and  k  *  .13,  but  to  a  lesser  extent.  In  these 
cases  the  calculated  lift  and  drag  characteristics.  Figures  19(b) 
and  (c),  respectively,  are  in  very  good  agreement  with  experi¬ 
ment.  The  tendency  for  the  measured  lift  to  peak  at  about  a  * 
30o  is  also  shown  in  the  calculated  results.  These  calculations 
were  not  continued  at  0^-^  for  a  sufficient  time  to  enable 
"steady-state"  conditions  to  be  reached.  As  stated  earlier, 
because  of  problems  with  the  amalgamation  routine,  a  strong 
vortex  core  did  not  "condense"  for  these  cases  and  so  the  upper 
surface  suction  peak  seen  in  the  experimental  measurements  (47) 
did  not  materialize;  rather,  a  smeared  suction  peak  appeared 
because  of  the  more  diffuse  region  of  shed  vorticity.  Conse¬ 
quently,  although  the  integrated  lift  and  drag  are  in  good  agree¬ 
ment,  the  pitching  moment  characteristic  (not  shown)  is  not 
satisfactory  at  this  time.  This  upper  surface  suction  peak, 
which  is  associated  with  a  reversed  flow  region  under  the  vortex, 
was,  in  fact,  computed  in  earlier  preliminary  calculations  (see 
Figure  13)  involving  a  higher  pitch  rate,  k  ■  .175.  In  this  case 
a  vortex  core  condensed  early  in  the  calculation  (Figure  13(a)). 

Overall  these  calculations  are  very  encouraging  and  with 
some  refinement  in  the  vortex  amalgamation  procedure,  it  is 
anticipated  that  the  details  of  the  unsteady  pressure  distribu¬ 
tions  will  be  achieved. 

Although  the  two-dimensional  pilot  program  was  generated 
primarily  as  a  tool  to  examine  the  behavior  of  various  parts  of 
the  dynamic  separation  calculation,  it  has  shown  considerable 
promise  as  a  general  purpose  code  for  two-dimensional  calcula¬ 
tions.  Earlier  examples  (e.g..  Figures  10,  11  and  12)  demon¬ 
strated  a  capability  to  compute  base  pressures  and  drag  coeffi¬ 
cients  of  blunt  sections  using  an  impulsive  start.  An  extension 
of  this  to  compute  spoiler  characteristics  has  also  been  briefly 
examined.  Figure  20  shows  computed  wake  configurations  at  two 
steps.  This  is  for  the  case  of  a  spoiler  deflected  30°  on  an 
airfoil  at  a  «  80.  The  final  base  pressure  and  integrated  lift. 
Figure  20  (b)  and  (c),  respectively,  are  in  good  agreement  with 
experimental  measurements  (48).  The  calculated  values  represent 
an  average  value  over  the  last  few  time  steps  as  the  solution  had 
started  to  oscillate.  The  amplitude  on  is  about  .1  but  the 
calculation  ought  to  be  continued  for  a  longer  time  to  examine 
whether  a  pattern  between  upper  and  lower  vortex  formation  is 
established.  This  application  could  be  extended  further  to 
examine  pitch  rates  and,  with  a  straightforward  extension  of  the 
code,  rates  of  spoiler  deployment.  Such  an  extension,  involving 
relative  motion  between  parts  of  the  configuration,  would  also 
allow  treatment  of  pitching  airfoils  between  channel  walls  to 
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assess  the  effects  of  unsteady  blockage  in  pitch-up  or  oscilla¬ 
tory  experiments. 

While  the  two-dimensional  program  has  been  used  to  examine 
and  develop  the  various  routines  required  for  the  coupled  dynamic 
separation  calculations,  the  three-dimensional  code  development 
has  been  following  closely  behind.  The  unsteady  boundary  layer 
calculation — which  is  performed  along  computed  surface  (external) 
streamlines  at  each  time  step — has  been  fully  coupled  with  the 
unsteady  inviscid  program.  Test  cases  have  been  performed  and 
compared  with  experimental  data  from  the  DFVLR-AVA  in  Gottingen. 
These  experiments  were  conducted  as  part  of  a  cooperative  agree¬ 
ment  between  the  DFVLR  Institute  of  Aeroelasticity/West  Germany 
and  NASA  Langley  Research  Center.  Figure  21  compares  the  calcu¬ 
lated  and  measured  real  and  imaginary  pressure  distributions  at  a 
70%  spanwise  station  on  an  AR  =  4  rectangular  wing  undergoing 
pitch  oscillation  about  the  quarter  chord  with  <*  *  7.9°  +  1.0 
sin  (.2t).  (Reynolds  number  is  1.35  x  10  .)  The  potential  flow 
solution  is  also  shown  to  indicate  the  extent  of  the  viscous 
correction.  The  complete  solution  is  in  very  good  agreement  with 
the  measurements.  This  is  still  true  for  the  condition,  a  «  12 
+  l.Oo  sin  (.3t),  which  is  approaching  the  condition  of  dynamic 
stall  onset;  a  pressure  deviation  is  apparent  near  the  leading 
edge.  Work  is  continuing  on  further  development  of  the  three- 
dimensional  method,  incorporating  the  techniques  that  are  being 
examined  in  the  two-dimensional  pilot  code. 
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Figure  21(a)  ,4i)  .  Comparison  of  Chordwise  Pressure  Distribution 

between  Computed  and  DFVLR  Data  (an  =  7.9°,  n 


Figure  21(a)  (jLi)  .  Comparison  of  Chordwise  Pressure  Distribution  (Imaginary  Part) 

at  y/s  =  0.70  between  Computed  and  DFVLR  Data  (a  =  7.9°,  a.  = 
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Figure  21 (b) . 


Comparison  of  Chordwise  Pressure  Distribution 
between  Computed  and  DFVLR  Data  («  =  12°, 

=  1.07;  k  =  .3,  Rn  =  1.35  x  10^;  °  1 

Spanwise  Location,  n  *  0.7. 


5.0  CONCLUSIONS 


A  system  of  routines  has  been  developed  to  couple  an 
unsteady  time-stepping  potential  flow  panel  method  with  an  exten¬ 
sive  separation  model  and  an  unsteady  boundary  layer  code.  The 
routines  include  treatment  of  the  growth  of  a  multiple  sheet 
dynamic  wake  model  and  also  the  movement  of  the  separation  loca¬ 
tion  with  time.  Preliminary  checkout  of  the  routines  using  a 
simplified  pilot  code  show  encouraging  results  for  conditions 
approaching  the  onset  of  dynamic  stall.  These  test  calculations 
have  shown  that  the  current  boundary  layer  procedure  is  capable 
of  predicting  the  unsteady  boundary  layer  development  and  it 
should  be  adequate  for  analyzing  dynamic  stall. 

Calculations  for  constant-rate  pitch-up  to  about  €0°  angle 
of  attack  have  shown  encouraging  agreement  with  experimental 
measurements  but  have  uncovered  a  weakness  in  the  vortex  amalga¬ 
mation  routine  in  the  pilot  program  Further  refinements  are 
therefore  planned  for  the  vortex  wake  treatment.  Work  is  con¬ 
tinuing  to  install  the  pilot  program  routines  in  the  three- 
dimensional  program  in  order  to  continue  the  investigation  of  the 
coupled  viscous/inviscid  approach  to  dynamic  separated  flow  cal¬ 
culation. 
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